set linesize 140
set rmsg on
set scheme s1mono
/******************************************************************************
	Project		:	Army Service in the All-Volunteer Era
	Author(s)	:	Kyle Greenberg	(kyle.greenberg@westpoint.edu)
					Matthew Gudgeon (matthew.gudgeon@westpoint.edu)
					Adam Isen 		(Adam.Isen@treasury.gov)
					Corbin Miller 	(Corbin.Miller@treasury.gov)
					Rich Patterson 	(rich_patterson@byu.edu)
	File Name	:	cumulative.do
	Description	:	This file estimates and outputs results for average outcomes
					over 0-19 and 11-19 years after application.
*******************************************************************************/

*---- if not running master, set up file structure ----*
*change directory to where programs and subfolders are stored
*cd ""

if "${raw}"==""		global raw		"raw/"
if "${data}"==""	global data		"data/"
if "${output}"==""	global output	"output/"

cap mkdir "${output}"
cap mkdir "${output}cumulative"

*---- select subprograms ----*
local avg_earn_emp_table	1
local compare_ols			1
local rf_wage_plot			1
local rf_ed_plot			1
local rf_emp_plot			1
local complier_destinies	1
local robust				1
local balanced_panel		1
local occupations			1

*----------------------------------------------------------------*
*---- average effect of enlisting on earnings and employment ----*
*----------------------------------------------------------------*
if `avg_earn_emp_table'==1{
	*initialize output file
	clear
	set obs 1
	foreach strata in all black white {
		foreach outcome in wages lnwages anywages anyf1098T any_ed school_q {
			foreach cut in 31 50 {
				foreach range in 0_19 11_19 {
					foreach x in b se m N {
						qui gen `x'_`cut'_`range'_`outcome'_`strata' = .
					}
				}
			}
		}
	}
	save "${output}cumulative/results_avg_earn_emp_ed", replace

	foreach strata in all black white {
		foreach outcome in wages lnwages anywages anyf1098T any_ed school_q {
			use "${data}cumulative-analysis" if `strata'==1, clear
			rename (everanyf1098T school_q) (any_ed_0_19_mean school_q_0_19_mean)
			gen any_ed_11_19_mean = any_ed_0_19_mean
			gen school_q_11_19_mean = school_q_0_19_mean

			*AFQT=31 cutoff, 0-19 yr cumulative estimates (weighted average specification) with all cohorts
			ivregress 2sls `outcome'_0_19_mean k31 k31_2 instk31 instk31_2 i.quarterFE `controls' (access=inst31) [aw=yrs_observed_0_19] if inrange(firstafqt,12,49) & inrange(firstyear,1990,2011), robust
			local b_31_0_19 = _b[access]
			local se_31_0_19 = _se[access]
			local N_31_0_19 = e(N)
			sum `outcome_reg' [aw=yrs_observed_0_19] if e(sample)
			local m_31_0_19 = r(mean)
		
			*AFQT=31 cutoff, 11-19 yr cumulative estimates (weighted average specification) with 1990-2007 cohorts
			ivregress 2sls `outcome'_11_19_mean k31 k31_2 instk31 instk31_2 i.quarterFE `controls' (access=inst31) [aw=yrs_observed_11_19] if inrange(firstafqt,12,49) & inrange(firstyear,1990,2007), robust
			local b_31_11_19 = _b[access]
			local se_31_11_19 = _se[access]
			local N_31_11_19 = e(N)
			sum `outcome'_11_19_mean [aw=yrs_observed_11_19] if e(sample)
			local m_31_11_19 = r(mean)
			
			*AFQT=50 cutoff, 0-19 yr cumulative estimates (weighted average specification) with all cohorts
			ivregress 2sls `outcome'_0_19_mean k50 k50_2 instk50 instk50_2 i.quarterFE `controls' (access=inst50) [aw=yrs_observed_0_19] if inrange(firstafqt,31,68) & inrange(firstyear,1990,2011), robust
			local b_50_0_19 = _b[access]
			local se_50_0_19 = _se[access]
			local N_50_0_19 = e(N)
			sum `outcome_reg' [aw=yrs_observed_0_19] if e(sample)
			local m_50_0_19 = r(mean)
		
			*AFQT=50 cutoff, 11-19 yr cumulative estimates (weighted average specification) with 1990-2007 cohorts
			ivregress 2sls `outcome'_11_19_mean k50 k50_2 instk50 instk50_2 i.quarterFE `controls' (access=inst50) [aw=yrs_observed_11_19] if inrange(firstafqt,31,68) & inrange(firstyear,1990,2007), robust
			local b_50_11_19 = _b[access]
			local se_50_11_19 = _se[access]
			local N_50_11_19 = e(N)
			sum `outcome'_11_19_mean [aw=yrs_observed_11_19] if e(sample)
			local m_50_11_19 = r(mean)
			
			use "${output}cumulative/results_avg_earn_emp_ed", clear
			foreach cut in 31 50{
				foreach range in 0_19 11_19 {
					foreach x in b se m N {
						di "`x'_`cut'_`range'_`outcome'_`strata' = ``x'_`cut'_`range''"
						replace `x'_`cut'_`range'_`outcome'_`strata' = ``x'_`cut'_`range''
					}
				}
			}
			save "${output}cumulative/results_avg_earn_emp_ed", replace
		}
	}

	*---- write results to table ----*
	use "${output}cumulative/results_avg_earn_emp_ed", clear
	
	*create significance stars
	foreach strata in all black white {
		foreach outcome in wages lnwages anywages anyf1098T any_ed school_q {
			foreach cut in 31 50{
				foreach range in 0_19 11_19 {
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') < 1.645	local star_`cut'_`range'_`outcome'_`strata' = "" 
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') >= 1.645	local star_`cut'_`range'_`outcome'_`strata' = "*" 
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') >= 1.96	local star_`cut'_`range'_`outcome'_`strata' = "**" 
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') >= 2.576	local star_`cut'_`range'_`outcome'_`strata' = "***" 
				}
			}
		}
	}

	file open c using "${output}cumulative/avg_earn_emp_all.tex", write replace
	file write c "\begin{tabular}{l*{4}c}" _n
	file write c "\hline\hline" _n
	file write c "& \multicolumn{2}{c}{31 AFQT Cutoff} & \multicolumn{2}{c}{50 AFQT Cutoff} \\" _n
	file write c "& 0-19 & 11-19  &  0-19 &11-19   \\" _n
	file write c "& Yrs Since &  Yrs Since &  Yrs Since &  Yrs Since   \\" _n
	file write c "& (1) & (2) & (3)  & (4)  \\ \hline" _n "\\" _n
	file write c "& \multicolumn{4}{c}{\underline{Panel (a): Earnings}} \\ \\" _n
	file write c "Enlist " "&" %5.0fc (b_31_0_19_wages_all) "`star_31_0_19_wages_all'" "&" %5.0fc (b_31_11_19_wages_all) "`star_31_11_19_wages_all'" "&" %5.0fc (b_50_0_19_wages_all) "`star_50_0_19_wages_all'" "&" %5.0fc (b_50_11_19_wages_all) "`star_50_11_19_wages_all'" "\\" _n
	file write c "& (" %5.0fc (se_31_0_19_wages_all) ") & ("  %5.0fc (se_31_11_19_wages_all) ") & ("  %5.0fc (se_50_0_19_wages_all) ") & ("  %5.0fc (se_50_11_19_wages_all) ") \\" _n
	file write c "Dep. Var Mean " "&" %6.0fc (m_31_0_19_wages_all) "&" %6.0fc (m_31_11_19_wages_all) "&" %6.0fc (m_50_0_19_wages_all) "&" %6.0fc (m_50_11_19_wages_all) "\\" _n "\\" _n
	file write c "& \multicolumn{4}{c}{\underline{Panel (b): Log Earnings}} \\ \\" _n
	file write c "Enlist " "&" %5.3f (b_31_0_19_lnwages_all) "`star_31_0_19_lnwages_all'" "&" %5.3f (b_31_11_19_lnwages_all) "`star_31_11_19_lnwages_all'" "&" %5.3f (b_50_0_19_lnwages_all) "`star_50_0_19_lnwages_all'" "&" %5.3f (b_50_11_19_lnwages_all) "`star_50_11_19_lnwages_all'" "\\" _n
	file write c "& (" %5.3f (se_31_0_19_lnwages_all) ") & ("  %5.3f (se_31_11_19_lnwages_all) ") & ("  %5.3f (se_50_0_19_lnwages_all) ") & ("  %5.3f (se_50_11_19_lnwages_all) ") \\" _n
	file write c "Dep. Var Mean " "&" %5.3f (m_31_0_19_lnwages_all) "&" %5.3f (m_31_11_19_lnwages_all) "&" %5.3f (m_50_0_19_lnwages_all) "&" %5.3f (m_50_11_19_lnwages_all) "\\" _n "\\" _n
	file write c "& \multicolumn{4}{c}{\underline{Panel (c): Employment}} \\ \\" _n
	file write c "Enlist &" %5.3f (b_31_0_19_anywages_all) "`star_31_0_19_anywages_all'" "&" %5.3f (b_31_11_19_anywages_all) "`star_31_11_19_anywages_all'" "&" %5.3f (b_50_0_19_anywages_all) "`star_50_0_19_anywages_all'" "&" %5.3f (b_50_11_19_anywages_all) "`star_50_11_19_anywages_all'" "\\" _n
	file write c "& (" %5.3f (se_31_0_19_anywages_all) ") & ("  %5.3f (se_31_11_19_anywages_all) ") & ("  %5.3f (se_50_0_19_anywages_all) ") & ("  %5.3f (se_50_11_19_anywages_all) ") \\" _n
	file write c "Dep. Var Mean " "&" %5.3f (m_31_0_19_anywages_all) "&" %5.3f (m_31_11_19_anywages_all) "&" %5.3f (m_50_0_19_anywages_all) "&" %5.3f (m_50_11_19_anywages_all) "\\" _n "\\" _n
	file write c "Observations " "&" %9.0fc (N_31_0_19_wages_all) "&" %9.0fc (N_31_11_19_wages_all) "&" %9.0fc (N_50_0_19_wages_all) "&" %9.0fc (N_50_11_19_wages_all) "\\" _n
	file write c "\hline \hline" _n "\end{tabular}" _n
	file close c
	
	*calculate p-values
	foreach outcome in wages lnwages anywages {
		foreach cut in 31 50{
			foreach range in 0_19 11_19 {
				local t_`cut'_`range'_`outcome'_bw = (b_`cut'_`range'_`outcome'_black - b_`cut'_`range'_`outcome'_white) / ((se_`cut'_`range'_`outcome'_black^2 + se_`cut'_`range'_`outcome'_white^2)^.5)
				local p_`cut'_`range'_`outcome'_bw = 2*ttail(N_`cut'_`range'_`outcome'_black+N_`cut'_`range'_`outcome'_white,abs(`t_`cut'_`range'_`outcome'_bw'))
			}
		}
	}
	
	file open c using "${output}cumulative/avg_earn_emp_race.tex", write replace
	file write c "\begin{tabular}{l*{4}c}" _n "\\ \hline\hline" _n
	file write c "& \multicolumn{2}{c}{\underline{0-19 Years}} & \multicolumn{2}{c}{\underline{11-19 Years}} \\" _n
	file write c "& Black &  White & Black &  White    \\ " _n
	file write c "& (1) & (2) & (3) & (4) \\ \cline{2-5} \\" _n
	file write c "& \multicolumn{4}{c}{\underline{Panel (a): Average Earnings, 31 Cutoff}} \\ \\" _n
	file write c "Enlist" "&" %5.0fc (b_31_0_19_wages_black) "`star_31_0_19_wages_black'" "&" %5.0fc (b_31_0_19_wages_white) "`star_31_0_19_wages_white'" "&" %5.0fc (b_31_11_19_wages_black) "`star_31_11_19_wages_black'" "&" %5.0fc (b_31_11_19_wages_black) "`star_31_11_19_wages_black'" "\\" _n
	file write c "& (" %5.0fc (se_31_0_19_wages_black) ") & (" %5.0fc (se_31_0_19_wages_white) ") & (" %5.0fc (se_31_11_19_wages_black) ") & (" %5.0fc (se_31_11_19_wages_white) ") \\" _n "\\" _n
	file write c "Observations " "&" %7.0fc (N_31_0_19_wages_black) "&" %7.0fc (N_31_0_19_wages_white) "&" %7.0fc (N_31_11_19_wages_black) "&" %7.0fc (N_31_11_19_wages_white) "\\" _n
	file write c "Dep. Var Mean " "&" %6.0fc (m_31_0_19_wages_black) "&" %6.0fc (m_31_0_19_wages_white) "&" %6.0fc (m_31_11_19_wages_black) "&" %6.0fc (m_31_11_19_wages_white) "\\" _n
	file write c "P-value for Equivalence & \multicolumn{2}{c}{" %5.3f (`p_31_0_19_wages_bw') "} & \multicolumn{2}{c}{" %5.3f (`p_31_11_19_wages_bw') "} \\" _n 
	file write c "\\" _n
	file write c "& \multicolumn{4}{c}{\underline{Panel (b): Average Earnings, 50 Cutoff}} \\ \\" _n
	file write c "Enlist" "&" %5.0fc (b_50_0_19_wages_black) "`star_50_0_19_wages_black'" "&" %5.0fc (b_50_0_19_wages_white) "`star_50_0_19_wages_white'" "&" %5.0fc (b_50_11_19_wages_black) "`star_50_11_19_wages_black'" "&" %5.0fc (b_50_11_19_wages_black) "`star_50_11_19_wages_black'" "\\" _n
	file write c "& (" %5.0fc (se_50_0_19_wages_black) ") & (" %5.0fc (se_50_0_19_wages_white) ") & (" %5.0fc (se_50_11_19_wages_black) ") & (" %5.0fc (se_50_11_19_wages_white) ") \\" _n "\\" _n
	file write c "Observations " "&" %7.0fc (N_50_0_19_wages_black) "&" %7.0fc (N_50_0_19_wages_white) "&" %7.0fc (N_50_11_19_wages_black) "&" %7.0fc (N_50_11_19_wages_white) "\\" _n
	file write c "Dep. Var Mean " "&" %6.0fc (m_50_0_19_wages_black) "&" %6.0fc (m_50_0_19_wages_white) "&" %6.0fc (m_50_11_19_wages_black) "&" %6.0fc (m_50_11_19_wages_white) "\\" _n
	file write c "P-value for Equivalence & \multicolumn{2}{c}{" %5.3f (`p_50_0_19_wages_bw') "} & \multicolumn{2}{c}{" %5.3f (`p_50_11_19_wages_bw') "} \\" _n 
	file write c "\\" _n
	file write c "& \multicolumn{4}{c}{\underline{Panel (c): Average Employment, 31 Cutoff}} \\ \\" _n
	file write c "Enlist" "&" %5.3f (b_31_0_19_anywages_black) "`star_31_0_19_anywages_black'" "&" %5.3f (b_31_0_19_anywages_white) "`star_31_0_19_anywages_white'" "&" %5.3f (b_31_11_19_anywages_black) "`star_31_11_19_anywages_black'" "&" %5.3f (b_31_11_19_anywages_black) "`star_31_11_19_anywages_black'" "\\" _n
	file write c "& (" %5.3f (se_31_0_19_anywages_black) ") & (" %5.3f (se_31_0_19_anywages_white) ") & (" %5.3f (se_31_11_19_anywages_black) ") & (" %5.3f (se_31_11_19_anywages_white) ") \\" _n "\\" _n
	file write c "Observations " "&" %7.0fc (N_31_0_19_anywages_black) "&" %7.0fc (N_31_0_19_anywages_white) "&" %7.0fc (N_31_11_19_anywages_black) "&" %7.0fc (N_31_11_19_anywages_white) "\\" _n
	file write c "Dep. Var Mean " "&" %5.3f (m_31_0_19_anywages_black) "&" %5.3f (m_31_0_19_anywages_white) "&" %5.3f (m_31_11_19_anywages_black) "&" %5.3f (m_31_11_19_anywages_white) "\\" _n
	file write c "P-value for Equivalence & \multicolumn{2}{c}{" %5.3f (`p_31_0_19_anywages_bw') "} & \multicolumn{2}{c}{" %5.3f (`p_31_11_19_anywages_bw') "} \\" _n 
	file write c "\\" _n
	file write c "& \multicolumn{4}{c}{\underline{Panel (d): Average Employment, 50 Cutoff}} \\ \\" _n
	file write c "Enlist" "&" %5.3f (b_50_0_19_anywages_black) "`star_50_0_19_anywages_black'" "&" %5.3f (b_50_0_19_anywages_white) "`star_50_0_19_anywages_white'" "&" %5.3f (b_50_11_19_anywages_black) "`star_50_11_19_anywages_black'" "&" %5.3f (b_50_11_19_anywages_black) "`star_50_11_19_anywages_black'" "\\" _n
	file write c "& (" %5.3f (se_50_0_19_anywages_black) ") & (" %5.3f (se_50_0_19_anywages_white) ") & (" %5.3f (se_50_11_19_anywages_black) ") & (" %5.3f (se_50_11_19_anywages_white) ") \\" _n "\\" _n
	file write c "Observations " "&" %7.0fc (N_50_0_19_anywages_black) "&" %7.0fc (N_50_0_19_anywages_white) "&" %7.0fc (N_50_11_19_anywages_black) "&" %7.0fc (N_50_11_19_anywages_white) "\\" _n
	file write c "Dep. Var Mean " "&" %5.3f (m_50_0_19_anywages_black) "&" %5.3f (m_50_0_19_anywages_white) "&" %5.3f (m_50_11_19_anywages_black) "&" %5.3f (m_50_11_19_anywages_white) "\\" _n
	file write c "P-value for Equivalence & \multicolumn{2}{c}{" %5.3f (`p_50_0_19_anywages_bw') "} & \multicolumn{2}{c}{" %5.3f (`p_50_11_19_anywages_bw') "} \\" _n 
	file write c "\hline \hline" _n "\end{tabular}" _n
	file close c

	*table a.14 panel a
	file open c "${output}/cumulative/table_a_14_panel_a.tex", write replace
	file write c "Enlist " 
	file write c "&" %5.0fc (b_31_0_19_school_q_black) `star_31_0_19_school_q_black' 
	file write c "&" %5.0fc (b_50_0_19_school_q_black) `star_50_0_19_school_q_black' 
	file write c "&" %5.0fc (b_31_0_19_school_q_white) `star_31_0_19_school_q_white' 
	file write c "&" %5.0fc (b_50_0_19_school_q_white) `star_50_0_19_school_q_white'
	file write c "\\" _n
	file write c " " 
	file write c "& (" %5.0fc (se_31_0_19_school_q_black)
	file write c ") & (" %5.0fc (se_50_0_19_school_q_black)
	file write c ") & (&" %5.0fc (se_31_0_19_school_q_white)
	file write c ") & (&" %5.0fc (se_50_0_19_school_q_white)
	file write c ") \\" _n
	file close c

}	//end if
else di "*---- skipping earnings and employment average effects ----*"

*---------------------*
*---- Compare OLS ----*
*---------------------*
if `compare_ols'==1{
	*initialize output file
	clear
	set obs 1
	foreach strata in all black white {
		foreach outcome in wages lnwages anywages {
			foreach cut in 31 50{
				foreach range in 0_19 11_19 {
					foreach x in b se m N {
						qui gen `x'_ols_`cut'_`range'_`outcome'_`strata' = .
					}
				}
			}
		}
	}
	save "${output}cumulative/results_ols_earn", replace

	use "${data}cumulative-analysis", clear
	
	*AFQT=31 cutoff, 0-19 yr cumulative estimates (weighted average specification) with all cohorts
	reg wages_0_19_mean access i.firstafqt i.quarterFE [aw=yrs_observed_0_19] if inrange(firstafqt,12,49) & inrange(firstyear,1990,2011), robust
	local b_31_0_19 = _b[access]
	local se_31_0_19 = _se[access]
	local N_31_0_19 = e(N)
	sum wages_0_19_mean [aw=yrs_observed_0_19] if e(sample)
	local m_31_0_19 = r(mean)

	*AFQT=31 cutoff, 11-19 yr cumulative estimates (weighted average specification) with 1990-2007 cohorts
	reg wages_11_19_mean access i.firstafqt i.quarterFE [aw=yrs_observed_11_19] if inrange(firstafqt,12,49) & inrange(firstyear,1990,2007), robust
	local b_31_11_19 = _b[access]
	local se_31_11_19 = _se[access]
	local N_31_11_19 = e(N)
	sum wages_11_19_mean [aw=yrs_observed_11_19] if e(sample)
	local m_31_11_19 = r(mean)
	
	*AFQT=50 cutoff, 0-19 yr cumulative estimates (weighted average specification) with all cohorts
	reg wages_0_19_mean access i.firstafqt i.quarterFE [aw=yrs_observed_0_19] if inrange(firstafqt,31,68) & inrange(firstyear,1990,2011), robust
	local b_50_0_19 = _b[access]
	local se_50_0_19 = _se[access]
	local N_50_0_19 = e(N)
	sum wages_0_19_mean [aw=yrs_observed_0_19] if e(sample)
	local m_50_0_19 = r(mean)

	*AFQT=50 cutoff, 11-19 yr cumulative estimates (weighted average specification) with 1990-2007 cohorts
	reg wages_11_19_mean access i.firstafqt i.quarterFE [aw=yrs_observed_11_19] if inrange(firstafqt,31,68) & inrange(firstyear,1990,2007), robust
	local b_50_11_19 = _b[access]
	local se_50_11_19 = _se[access]
	local N_50_11_19 = e(N)
	sum wages_11_19_mean [aw=yrs_observed_11_19] if e(sample)
	local m_50_11_19 = r(mean)
	
	use "${output}cumulative/results_ols_earn", clear
	foreach cut in 31 50 {
		foreach range in 0_19 11_19 {
			foreach x in b se m N {
				di "`x'_`cut'_`range'_wages_all = ``x'_`cut'_`range''"
				replace `x'_ols_`cut'_`range'_wages_all = ``x'_`cut'_`range''
			}
		}
	}
	save "${output}cumulative/results_ols_earn", replace

	file open c using "${output}cumulative/compare_ols.tex", write replace
	file write c "\begin{tabular}{l*{4}c}" _n "\hline\hline" _n
	file write c "& \multicolumn{2}{c}{31 AFQT Cutoff} &  \multicolumn{2}{c}{50 AFQT Cutoff} \\" _n
	file write c "& 0-19 & 11-19  & 0-19 & 11-19  \\" _n
	file write c "& Yrs Since &  Yrs Since & Yrs Since&  Yrs Since\\" _n
	file write c "& (1) & (2) &  (3)  & (4)  \\ \hline" _n
	file write c "\\" _n "& \multicolumn{4}{c}{\underline{Panel (a): Primary (2SLS RD) Estimates}} \\ \\" _n
	use "${output}cumulative/results_avg_earn_emp_ed", clear
	*create significance stars
	foreach strata in all {
		foreach outcome in wages {
			foreach cut in 31 50{
				foreach range in 0_19 11_19 {
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') < 1.645	local star_`cut'_`range'_`outcome'_`strata' = "" 
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') >= 1.645	local star_`cut'_`range'_`outcome'_`strata' = "*" 
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') >= 1.96	local star_`cut'_`range'_`outcome'_`strata' = "**" 
					if abs(b_`cut'_`range'_`outcome'_`strata'/se_`cut'_`range'_`outcome'_`strata') >= 2.576	local star_`cut'_`range'_`outcome'_`strata' = "***" 
				}
			}
		}
	}
	file write c "Enlist " "&" %5.0fc (b_31_0_19_wages_all) "`star_31_0_19_wages_all'" "&" %5.0fc  (b_31_11_19_wages_all) "`star_31_11_19_wages_all'" "&" %5.0fc (b_50_0_19_wages_all) "`star_50_0_19_wages_all'" "&" %5.0fc (b_50_11_19_wages_all) "`star_50_11_19_wages_all'" "\\" _n
	file write c "" "& (" %5.0fc (se_31_0_19_wages_all) ") & (" %5.0fc (se_31_11_19_wages_all) ") & (" %5.0fc (se_50_0_19_wages_all) ") & (" %5.0fc (se_50_11_19_wages_all) ") \\" _n
	file write c "\\" _n "& \multicolumn{4}{c}{\underline{Panel (b): OLS Estimates}} \\ \\" _n
	use "${output}cumulative/results_ols_earn", clear
	foreach cut in 31 50{
		foreach range in 0_19 11_19{
			if abs(b_ols_`cut'_`range'_wages_all/se_ols_`cut'_`range'_wages_all) < 1.645	local star_`cut'_`range'_wages_all = "" 
			if abs(b_ols_`cut'_`range'_wages_all/se_ols_`cut'_`range'_wages_all) >= 1.645	local star_`cut'_`range'_wages_all = "*" 
			if abs(b_ols_`cut'_`range'_wages_all/se_ols_`cut'_`range'_wages_all) >= 1.96	local star_`cut'_`range'_wages_all = "**" 
			if abs(b_ols_`cut'_`range'_wages_all/se_ols_`cut'_`range'_wages_all) >= 2.576	local star_`cut'_`range'_wages_all = "***" 
		}
	}
	file write c "Enlist " "&" %5.0fc (b_ols_31_0_19_wages_all) "`star_31_0_19_wages_all'" "&" %5.0fc (b_ols_31_11_19_wages_all) "`star_31_11_19_wages_all'" "&" %5.0fc (b_ols_50_0_19_wages_all) "`star_50_0_19_wages_all'" "&" %5.0fc (b_ols_50_11_19_wages_all) "`star_50_11_19_wages_all'" "\\" _n
	file write c "" "& (" %2.0fc (se_ols_31_0_19_wages_all) ") & (" %2.0fc (se_ols_31_11_19_wages_all) ") & (" %2.0fc (se_ols_31_0_19_wages_all) ") & (" %2.0fc (se_ols_50_11_19_wages_all) ") \\" _n
	file write c "\\ \hline" _n
	file write c "Dep. Var Mean " "&" %6.0fc (m_ols_31_0_19_wages_all) "&" %6.0fc (m_ols_31_11_19_wages_all) "&" %6.0fc (m_ols_50_0_19_wages_all) "&" %6.0fc (m_ols_50_11_19_wages_all) "\\" _n
	file write c "Observations " "&" %9.0fc (N_ols_31_0_19_wages_all) "&" %9.0fc (N_ols_31_11_19_wages_all) "&" %9.0fc (N_ols_50_0_19_wages_all) "&" %9.0fc (N_ols_50_11_19_wages_all) "\\" _n
	file write c "\hline \hline" _n "\end{tabular}" _n
	file close c

}	//end if
else di "*---- skipping OLS comparison ----*"

*-----------------------------* 
*---- reduced-form graphs ----*
*-----------------------------*

*---- Earnings ----*
if `rf_wage_plot'==1{
	foreach strata in all black white hisp {
		foreach range in 11_19 0_19 {  
			use "${data}cumulative-analysis" if `strata' == 1 & yrs_observed_`range' > 0 & !missing(yrs_observed_`range',wages_`range'_mean), clear
			replace firstafqt_fy=2004.75 if firstafqt_fy==3000
			*omit first quarter of FY2000 as base mean
			reg wages_`range'_mean ib40.quarterFE [aw=yrs_observed_`range'] if inrange(firstafqt,12,68)
			predict r if inrange(firstafqt,12,68), resid
			sum wages_`range'_mean [aw=yrs_observed_`range'] if e(sample)
			replace wages_`range'_mean = r + r(mean) if inrange(firstafqt,12,68)
			collapse (mean) wages_`range'_mean [aw=yrs_observed_`range'], by(firstafqt)
			tw (qfitci wages_`range'_mean firstafqt if inrange(firstafqt,12,30), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(qfitci wages_`range'_mean firstafqt if inrange(firstafqt,31,49), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(qfitci wages_`range'_mean firstafqt if inrange(firstafqt,50,68), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(scatter wages_`range'_mean firstafqt if inrange(firstafqt,12,68), msiz(medsmall) mc(black) mfc(black) mlw(vvthin)) ///
				,  legend(off) xline(30.5 49.5, lcolor(black) lpattern(dot)) ///
				ytitle("Earnings (Dollars)", size(small)) ylabel(, angle(0) format(%9.0fc))  ///
				xtitle("First AFQT Score on File") xlabel(15(5)65) xmtick(12(1)68)
			graph export "${output}cumulative/wages_rf_wgt_`strata'_`range'.pdf", as(pdf) replace
		}
	}
}	//end if
else di "*---- skipping reduced form wage plots ----*"

*---- Education ----*
if `rf_ed_plot'==1{
	foreach strata in all black white hisp {
		foreach range in 11_19 0_19 {  
			use "${data}cumulative-analysis" if `strata' == 1 & yrs_observed_`range' > 0 & !missing(yrs_observed_`range',wages_`range'_mean), clear
			replace firstafqt_fy=2004.75 if firstafqt_fy==3000
			reg anyf1098T_`range'_mean ib40.quarterFE [aw=yrs_observed_`range'] if inrange(firstafqt,12,68)
			predict r if inrange(firstafqt,12,68), resid
			sum anyf1098T_`range'_mean [aw=yrs_observed_`range'] if e(sample)
			replace anyf1098T_`range'_mean = r + r(mean) if inrange(firstafqt,12,68)
			collapse (mean) anyf1098T_`range'_mean [aw=yrs_observed_`range'], by(firstafqt)
			tw (qfitci anyf1098T_`range'_mean firstafqt if inrange(firstafqt,12,30), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(qfitci anyf1098T_`range'_mean firstafqt if inrange(firstafqt,31,49), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(qfitci anyf1098T_`range'_mean firstafqt if inrange(firstafqt,50,68), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(scatter anyf1098T_`range'_mean firstafqt if inrange(firstafqt,12,68), msiz(medsmall) mc(black) mfc(black) mlw(vvthin)) ///
				,  legend(off) xline(30.5 49.5, lcolor(black) lpattern(dot)) ///
				ytitle("Receives Form 1098-T", size(small)) ylabel(, angle(0) format(%5.2f))  ///
				xtitle("First AFQT Score on File") xlabel(15(5)65) xmtick(12(1)68)
			graph export "${output}cumulative/anyf1098T_rf_wgt_`strata'_`range'.pdf", as(pdf) replace
		}
	}
}	//end if
else di "*---- skipping reduced form education plots ----*"

*---- Employment ----*
if `rf_emp_plot'==1{
	foreach strata in all black white hisp {
		foreach range in 11_19 0_19 {  
			use "${data}cumulative-analysis" if `strata' == 1 & yrs_observed_`range' > 0 & !missing(yrs_observed_`range',wages_`range'_mean), clear
			replace firstafqt_fy=2004.75 if firstafqt_fy==3000
			reg anywages_`range'_mean ib40.quarterFE [aw=yrs_observed_`range'] if inrange(firstafqt,12,68)
			predict r if inrange(firstafqt,12,68), resid
			sum anywages_`range'_mean [aw=yrs_observed_`range'] if e(sample)
			replace anywages_`range'_mean = r + r(mean) if inrange(firstafqt,12,68)
			collapse (mean) anywages_`range'_mean [aw=yrs_observed_`range'], by(firstafqt)
			tw (qfitci anywages_`range'_mean firstafqt if inrange(firstafqt,12,30), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(qfitci anywages_`range'_mean firstafqt if inrange(firstafqt,31,49), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(qfitci anywages_`range'_mean firstafqt if inrange(firstafqt,50,68), clc(black) clw(medthick) fc(gs13) alc(white)) ///
				(scatter anywages_`range'_mean firstafqt if inrange(firstafqt,12,68), msiz(medsmall) mc(black) mfc(black) mlw(vvthin)) ///
				,  legend(off) xline(30.5 49.5, lcolor(black) lpattern(dot)) ///
				ytitle("Employed", size(small)) ylabel(, angle(0) format(%5.2f))  ///
				xtitle("First AFQT Score on File") xlabel(15(5)65) xmtick(12(1)68)
			graph export "${output}cumulative/anywages_rf_wgt_`strata'_`range'.pdf", as(pdf) replace
		}
	}
}	//end if
else di "*---- skipping reduced form employment plots ----*"

*---- complier destinies ----*
if `complier_destinies'==1{
	use "${data}cumulative-analysis", clear
	
	merge 1:1 pid using "${raw}mepcom-treasury", keep(1 3) nogen keepusing(access_service_1 access* active otherserv reserve min_access_service min_access_dt)

	gen oth_active = inlist(min_access_serv,"CR","FR","MR","NR")
	gen oth_reserve = inlist(min_access_serv,"CV","FG","FV","MV","NV")
	gen army_reserve = inlist(min_access_serv,"AG","AV")
	gen active_navy = min_access_serv=="NR" 
	gen active_af = min_access_serv=="FR" 
	gen active_cg = min_access_serv=="CR" 
	gen active_mc = min_access_serv=="MR" 
	gen army_reserve_guard = accessarmy_reserves | accessarmy_guard
	gen active_other = active & !accessarmy_active
	gen no_access = !access

	local access_list "no_access access active active_other reserve accessarmy_active active_mc active_navy active_af active_cg army_reserve_guard accessarmy_guard accessarmy_reserves oth_reserve"

	foreach strata in all black white{
		*loop through accession option (rows of the table). Save sample mean, analysis sample mean, and first stage regression in matricies.
		foreach x of local access_list{
			qui summ `x' if `strata'
			matrix `x'_mean = r(mean)
			qui summ `x' if inrange(firstafqt,12,68) & inrange(firstyear,1990,2011) & `strata'
			matrix `x'_anal_mean = r(mean)
			qui summ `x' if inrange(firstafqt,12,49) & inrange(firstyear,1990,2011) & `strata'
			matrix `x'_31_mean = r(mean)
			qui summ `x' if inrange(firstafqt,31,68) & inrange(firstyear,1990,2011) & `strata'
			matrix `x'_50_mean = r(mean)	
			reghdfe `x' inst31 k31 k31_2 instk31 instk31_2 if inrange(firstafqt,12,49) & inrange(firstyear,1990,2011) & `strata', a(quarterFE)
			matrix `x'_fs31 = r(table)
			reghdfe `x' inst50 k50 k50_2 instk50 instk50_2 if inrange(firstafqt,31,68) & inrange(firstyear,1990,2011) & `strata', a(quarterFE)
			matrix `x'_fs50 = r(table)
			
			local `x'_star_31 
			if `x'_fs31[4,1] < .1 local `x'_star_31 = "*"
			if `x'_fs31[4,1] < .05 local `x'_star_31 = "**"
			if `x'_fs31[4,1] < .01 local `x'_star_31 = "***"
			local `x'_star_50 
			if `x'_fs50[4,1] < .1 local `x'_star_50 = "*"
			if `x'_fs50[4,1] < .05 local `x'_star_50 = "**"
			if `x'_fs50[4,1] < .01 local `x'_star_50 = "***"	
		}

		*Save number of soldiers in data, analysis sample, 31 cutoff sample, and 50 cutoff sample.
		summ access if `strata'
		matrix all_N = r(N)
		summ access if inrange(firstyear,1990,2011) & inrange(firstafqt,12,68) & `strata'
		matrix analysis_N = r(N)
		summ access if inrange(firstyear,1990,2011) & inrange(firstafqt,12,49) & `strata'
		matrix fs31_N = r(N)
		summ access if inrange(firstyear,1990,2011) & inrange(firstafqt,31,68) & `strata'
		matrix fs50_N = r(N)

		foreach x of local access_list{
			foreach y in mean anal_mean 31_mean 50_mean fs31 fs50 {
				local `x'_`y' = `x'_`y'[1,1]
			}
			local `x'_star_fs31 
			if `x'_fs31[4,1] < .1 local `x'_star_fs31 = "*"
			if `x'_fs31[4,1] < .05 local `x'_star_fs31 = "**"
			if `x'_fs31[4,1] < .01 local `x'_star_fs31 = "***"
			local `x'_star_fs50 
			if `x'_fs50[4,1] < .1 local `x'_star_fs50 = "*"
			if `x'_fs50[4,1] < .05 local `x'_star_fs50 = "**"
			if `x'_fs50[4,1] < .01 local `x'_star_fs50 = "***"
		}

		file open x using "${output}cumulative/complier_destinies_`strata'.tex", write replace
		file write x "\multicolumn{1}{l}{Any Military}"_char(38) %5.4f (access_mean[1,1]) _char(38) %5.4f (access_31_mean[1,1]) _char(38) %5.4f (access_fs31[1,1]) ("`access_star_fs31'") _char(38) _char(38) %5.4f (access_50_mean[1,1]) _char(38) %5.4f (access_fs50[1,1]) ("`access_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (access_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (access_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "Any Active Duty"_char(38) %5.4f (active_mean[1,1]) _char(38) %5.4f (active_31_mean[1,1]) _char(38) %5.4f (active_fs31[1,1]) ("`active_star_fs31'") _char(38) _char(38) %5.4f (active_50_mean[1,1]) _char(38) %5.4f (active_fs50[1,1]) ("`active_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (active_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (active_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Army}" _char(38) %5.4f (accessarmy_active_mean[1,1]) _char(38) %5.4f (accessarmy_active_31_mean[1,1]) _char(38) %5.4f (accessarmy_active_fs31[1,1]) ("`accessarmy_active_star_fs31'") _char(38) _char(38) %5.4f (accessarmy_active_50_mean[1,1]) _char(38) %5.4f (accessarmy_active_fs50[1,1]) ("`accessarmy_active_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (accessarmy_active_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (accessarmy_active_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Marines}"_char(38) %5.4f (active_mc_mean[1,1]) _char(38) %5.4f (active_mc_31_mean[1,1]) _char(38) %5.4f (active_mc_fs31[1,1]) ("`active_mc_star_fs31'") _char(38) _char(38) %5.4f (active_mc_50_mean[1,1]) _char(38) %5.4f (active_mc_fs50[1,1]) ("`active_mc_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (active_mc_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (active_mc_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Navy}"_char(38) %5.4f (active_navy_mean[1,1]) _char(38) %5.4f (active_navy_31_mean[1,1]) _char(38) %5.4f (active_navy_fs31[1,1]) ("`active_navy_star_fs31'") _char(38) _char(38) %5.4f (active_navy_50_mean[1,1]) _char(38) %5.4f (active_navy_fs50[1,1]) ("`active_navy_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (active_navy_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (active_navy_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Air Force}"_char(38) %5.4f (active_af_mean[1,1]) _char(38) %5.4f (active_af_31_mean[1,1]) _char(38) %5.4f (active_af_fs31[1,1]) ("`active_af_star_fs31'") _char(38) _char(38) %5.4f (active_af_50_mean[1,1]) _char(38) %5.4f (active_af_fs50[1,1]) ("`active_af_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (active_af_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (active_af_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Coast Guard}"_char(38) %5.4f (active_cg_mean[1,1]) _char(38) %5.4f (active_cg_31_mean[1,1]) _char(38) %5.4f (active_cg_fs31[1,1]) ("`active_cg_star_fs31'") _char(38) _char(38) %5.4f (active_cg_50_mean[1,1]) _char(38) %5.4f (active_cg_fs50[1,1]) ("`active_cg_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (active_cg_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (active_cg_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Any Reserve Component}"_char(38) %5.4f (reserve_mean[1,1]) _char(38) %5.4f (reserve_31_mean[1,1]) _char(38) %5.4f (reserve_fs31[1,1]) ("`reserve_star_fs31'") _char(38) _char(38) %5.4f (reserve_50_mean[1,1]) _char(38) %5.4f (reserve_fs50[1,1]) ("`reserve_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (reserve_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (reserve_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Army Reserves}"_char(38) %5.4f (accessarmy_reserves_mean[1,1]) _char(38) %5.4f (accessarmy_reserves_31_mean[1,1]) _char(38) %5.3f (accessarmy_reserves_fs31[1,1]) ("`accessarmy_reserves_star_fs31'") _char(38) _char(38) %5.4f (accessarmy_reserves_50_mean[1,1]) _char(38) %5.4f (accessarmy_reserves_fs50[1,1]) ("`accessarmy_reserves_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (accessarmy_reserves_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (accessarmy_reserves_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Army National Guard}"_char(38) %5.4f (accessarmy_guard_mean[1,1]) _char(38) %5.4f (accessarmy_guard_31_mean[1,1]) _char(38) %5.4f (accessarmy_guard_fs31[1,1]) ("`accessarmy_guard_star_fs31'") _char(38) _char(38) %5.4f (accessarmy_guard_50_mean[1,1]) _char(38) %5.4f (accessarmy_guard_fs50[1,1]) ("`accessarmy_guard_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (accessarmy_guard_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (accessarmy_guard_fs50[2,1]) _char(41) _char(92) _char(92) _n
		
		file write x "\multicolumn{1}{l}{Other Reserve/Guard}"_char(38) %5.4f (oth_reserve_mean[1,1]) _char(38) %5.4f (oth_reserve_31_mean[1,1]) _char(38) %5.4f (oth_reserve_fs31[1,1]) ("`oth_reserve_star_fs31'") _char(38) _char(38) %5.4f (oth_reserve_50_mean[1,1]) _char(38) %5.4f (oth_reserve_fs50[1,1]) ("`oth_reserve_star_fs50'") _char(92) _char(92) _n
		file write x "" _char(38) _char(38) _char(38) %5.4f _char(40) (oth_reserve_fs31[2,1]) _char(41) _char(38) _char(38) _char(38) %5.4f _char(40) (oth_reserve_fs50[2,1]) _char(41) _char(92) _char(92) _n
		file write x "\hline" _char(92) _char(92)
		file write x _char(92) _char(92) "N" _char(38) (all_N[1,1]) _char(38) (fs31_N[1,1]) _char(38) _char(38) _char(38) (fs50_N[1,1]) _char(38) _char(92) _char(92) _n
		file write x "\hline\hline"
		file close x
	}
}	//end if
else di "*---- skipping complier destinies ----*"

*------------------------------------------------------------------------------------*
*---- Average Effects on Earnings 11-19 Years Post Application, Inference Checks ----*
*------------------------------------------------------------------------------------*
if `robust'==1{
	use "${data}cumulative-analysis", clear
	foreach n in 31 50 {
		gen tri`n' = 1-abs((firstafqt - `n'+.5)/19)
		gen tri`n'_cum_0_19=tri`n'* yrs_observed_0_19
		gen tri`n'_cum_11_19=tri`n'*yrs_observed_11_19
	}
	local outcome wages
	foreach strata in all black white{
		eststo clear
		foreach cut in 31 50{
			preserve
			if `cut'==31 keep if inrange(firstafqt,12,49) 
			if `cut'==50 keep if inrange(firstafqt,31,68)

			keep `strata' `outcome'_11_19_mean tri`cut'_cum_11_19 firstafqt k`cut' k`cut'_2 instk`cut' instk`cut'_2 quarterFE access inst`cut'  yrs_observed_11_19  firstyear
			sum tri`cut'_cum_11_19 if inrange(firstyear,1990,2007)
			keep if `strata' == 1 

			*2SLS
			eststo out_11_19_`cut': ivregress 2sls `outcome'_11_19_mean k`cut' k`cut'_2 instk`cut' instk`cut'_2 i.quarterFE (access=inst`cut') [aweight = yrs_observed_11_19] if inrange(firstyear,1990,2007) , robust
			qui summ `outcome'_11_19_mean if e(sample)
			estadd scalar depmean = r(mean)

			eststo out_11_19_c_`cut': ivregress 2sls `outcome'_11_19_mean k`cut' k`cut'_2 instk`cut' instk`cut'_2 i.quarterFE (access=inst`cut') [aweight = yrs_observed_11_19] if inrange(firstyear,1990,2007) , cluster(firstafqt) robust
			qui summ `outcome'_11_19_mean if e(sample)
			estadd scalar depmean = r(mean)

			drop access
			g access = inst`cut'
			*REDUCED FORMS
			eststo out_11_19_rf_`cut': regress `outcome'_11_19_mean access k`cut' k`cut'_2 instk`cut' instk`cut'_2 i.quarterFE [aweight = yrs_observed_11_19] if inrange(firstyear,1990,2007) , robust
			qui summ `outcome'_11_19_mean if e(sample)
			estadd scalar depmean = r(mean)

			eststo out_11_19_rfc_`cut': regress `outcome'_11_19_mean access k`cut' k`cut'_2 instk`cut' instk`cut'_2 i.quarterFE [aweight = yrs_observed_11_19] if inrange(firstyear,1990,2007) , robust cluster(firstafqt) 
			qui summ `outcome'_11_19_mean if e(sample)
			estadd scalar depmean = r(mean)

			*GROUPED MEANS
			matrix groupdata=J(68,18,.)
			local colg=1
			local rowg=1
			qui sum firstafqt if `strata'==1
			forvalues k=`r(min)'(1)`r(max)'{
				matrix groupdata[`k',`colg']=`k'
				*Include number of applicant-years with each AFQT score (so we can weight properly)
				qui sum yrs_observed_11_19 if firstafqt==`k' & `strata'==1
				matrix groupdata[`k',`colg'+1]=`r(sum)'
			}
			*Construct AFQT dummies
			qui sum firstafqt if `strata'==1
			forvalues k=`r(min)'(1)`r(max)'{
				gen afqt_g`k'=(firstafqt==`k')
			}
			*Construct covariate adjusted group effects
			local colg=3
			* Full Window, Quadratic, No Kernel
			qui xi: regress `outcome'_11_19_mean afqt_g* i.quarterFE [aweight = yrs_observed_11_19] if inrange(firstyear,1990,2007) & `strata', noconstant
			qui sum firstafqt
			forvalues k=`r(min)'(1)`r(max)'{
				matrix groupdata[`k',`colg'+2]=_b[afqt_g`k']
			}
			local colg=`colg'+4
			
			clear 
			qui svmat groupdata
										   
			rename groupdata1 firstafqt
			rename groupdata2 afqtn
			rename groupdata3 groupeffect_full_lin
			rename groupdata4 groupeffect_half_lin
			rename groupdata5 groupeffect_full_quad

			foreach n in 31 50 {
				gen inst`n'=(firstafqt>=`n')
				gen k`n'=firstafqt-`n'+.5
				gen instk`n'=inst`n'*k`n'      
				gen tri`n' = 1-abs((firstafqt - `n'+.5)/19)
				gen tri`n'_half = 1-abs((firstafqt - `n'+.5)/9)
				gen wt`n'=tri`n'*afqtn
				gen wt`n'_half=tri`n'_half*afqtn
				gen k`n'_2=(k`n')^2
				gen instk`n'_2=(instk`n')^2
			}

			*Group Means (Full Window, Quadratic, No Kernel)
			gen access = inst`cut'
			eststo out_11_19_gm_`cut': reg groupeffect_full_quad access k`cut' k`cut'_2 instk`cut' instk`cut'_2 [aweight = afqtn] if wt`cut'>0
			qui summ groupeffect_full_quad if e(strata)
			estadd scalar depmean = r(mean)

			restore
		}	// end cut loop
		
		noi esttab  out_11_19_31 out_11_19_c_31 out_11_19_rf_31 out_11_19_rfc_31 out_11_19_gm_31 out_11_19_50 out_11_19_c_50 out_11_19_rf_50 out_11_19_rfc_50 out_11_19_gm_50 using "${output}/cumulative/robust_earnings_`strata'.tex",	///
			cells(b(fmt(%9.0fc) star) se(par fmt(%9.0fc)) p(par([ ]) fmt(3)))  ///
			replace noobs nomtitle nodepvar nonumber plain fragment label ///	
			substitute(_ \_) style(tex) ///	
			mlabels(, none) collabels(, none) 	///			
			keep(access inst31 inst50 , relax) ///
			varlabels(access "Estimated Coef.") ///
			starlevels(* 0.1 ** 0.05 *** 0.01) varwidth(30) ///
			stats(N depmean, labels("\hline \\ Observations" "Dep. Var Mean") fmt(%16.0gc %9.0fc))
	}	//end strata loop
}	//end if

*-------------------------------------------------------------------------------------------------------------*
*---- average effects on earnings and employment 11-19 years after application: 90-99 application cohorts ----*
*-------------------------------------------------------------------------------------------------------------*
if `balanced_panel'==1{
	use "${data}cumulative-analysis", clear
	foreach strata in all black white {
		foreach outcome in wages anywages {
		eststo clear
			foreach cut in 31 50 {
				preserve
					if `cut'==31{
						keep if inrange(firstafqt,12,49) & `strata' == 1 
						keep `outcome'_0_19_mean `outcome'_11_19_mean `outcome'_15_19_mean k`cut' k`cut'_2 instk`cut' instk`cut'_2 quarterFE access inst`cut' yrs_observed_0_19 yrs_observed_11_19 yrs_observed_15_19 firstyear
					}
					if `cut'==50{
						keep if inrange(firstafqt,31,68) & `strata' == 1 
						keep `outcome'_0_19_mean `outcome'_11_19_mean `outcome'_15_19_mean k`cut' k`cut'_2 instk`cut' instk`cut'_2 quarterFE access inst`cut' yrs_observed_0_19 yrs_observed_11_19 yrs_observed_15_19 firstyear
					}
					
					eststo `outcome'_11_19_`cut': ivregress 2sls `outcome'_11_19_mean k`cut' k`cut'_2 instk`cut' instk`cut'_2 i.quarterFE (access=inst`cut') [aweight = yrs_observed_11_19] if inrange(firstyear,1990,1999), robust
					qui summ `outcome'_11_19_mean if e(sample)
					estadd scalar depmean = r(mean)
				restore
			}

			if inlist("`outcome'","wages","dis_tot") local formt %16.0gc
			else local formt 3

			noi esttab  `outcome'_11_19_31 `outcome'_11_19_50 ///
				using "${output}cumulative/balanced_panel_`outcome'_`strata'.tex",	///
				cells(b(fmt(`formt') star) se(par fmt(`formt')))  ///
				replace noobs nomtitle nodepvar nonumber plain fragment label ///	
				substitute(_ \_) style(tex) ///	
				mlabels(, none) collabels(, none) 	///			
				keep(access, relax) ///
				varlabels(access "Enlist (2SLS)") ///
				starlevels(* 0.1 ** 0.05 *** 0.01) varwidth(30) ///
				stats(N depmean, labels("\hline \\ Observations" "Dep. Var Mean") fmt(%16.0gc `formt'))
		}
	}
}	//end if
else di "*---- skipping robustness ----*"

*---------------------------------------------------------*
*---- occupations, deployment, combat injuries/deaths ----*
*---------------------------------------------------------*
if `occupations'==1{

	*---- clean dynamic wia ----*
	use "${raw}cas_mepcom_with_id_v1_nopii", clear
	keep if dodgrp == "HIW"
	g calendar_year = year(incidentdt_s)
	sort pid calendar_year 
	keep pid calendar_year 
	duplicates drop
	tsset pid calendar_year 
	g everwia = 1 
	tsfill, full
	bys pid (calendar_year): replace everwia = 1 if everwia[_n-1] == 1
	replace everwia = 0 if everwia == .
	tsset pid calendar_year
	tab calendar_year
	tempfile everwia
	save `everwia', replace

	*---- clean dynamic kia ----*
	use "${raw}cas_mepcom_with_id_v1_nopii", clear
	keep if dodgrp =="HDK" | dodgrp=="HDR" | dodgrp=="HDC" | dodgrp=="HDM"
	g calendar_year = year(incidentdt_s)
	sort pid calendar_year 
	keep pid calendar_year 
	duplicates drop
	tsset pid calendar_year 
	g everkia = 1 
	tsfill, full
	bys pid (calendar_year): replace everkia = 1 if everkia[_n-1] == 1
	replace everkia = 0 if everkia == .
	tsset pid calendar_year
	tab calendar_year
	tempfile everkia
	save `everkia', replace

	*---- clean up dynamic any hfp ----*
	use "${raw}pay-panel-cy-1992-2019_v2_v2", clear
	g calendar_year = year(ap_filedt_s)
	g everhfp = ap_y_hfp > 0 & ap_y_hfp < .
	keep pid calendar_year everhfp
	gsort pid calendar_year -everhfp
	by pid calendar_year: keep if _n == 1
	tsset pid calendar_year 
	tsfill, full
	bys pid (calendar_year): replace everhfp = 1 if everhfp[_n-1] == 1
	replace everhfp = 0 if everhfp == .
	tsset pid calendar_year 
	tab calendar_year
	tempfile everhfp
	save `everhfp', replace

	*---- create combat data ----*
	use "${data}cumulative-analysis", clear
	merge 1:1 pid using "${raw}mepcom-treasury", nogen keep(1 3) keepusing(init_pmos access_dt_s_1 firstafqtdt)
	gen month = month(access_dt_s_1)
	gen year = year(access_dt_s_1)
	gen mos = init_pmos
	replace mos = lower(mos)

	replace month = 1 if year< 1997
	replace year = 1997 if year < 1997
	replace month = 12 if year > 2010 & year < .
	replace year = 2010 if year > 2010 & year < .

	merge m:1 mos year month using "${raw}pmos-to-conspmos-1997-2010a-forsrb"
	drop month year
	tab _merge if accessarmy_active == 1
	keep if _merge == 1 | _merge == 3
	drop _merge

	gen combat_tc =(inlist(comcmf,11,13,14,15,18,19))
	replace combat_tc =1 if conspmos=="12b" 
	gen noncombat_tc = combat_tc != 1 
	replace noncombat_tc = 0 if conspmos == ""

	keep pid combat_tc noncombat_tc firstyear comcmf conspmos black white firstafqt accessarmy_active accessarmy inst31 k31 k31_2 instk31 instk31_2 inst50 k50 k50_2 instk50 instk50_2 quarterFE

	*---- merge in dynamic hfp ----*
	forval y=-1/19 {
		local ytitle = `y'
		if `y' == -1 local ytitle m1
		gen calendar_year = firstyear+`y'		
		merge m:1 pid calendar_year using `everhfp', nogen keep(1 3)
		merge m:1 pid calendar_year using `everwia', nogen keep(1 3)
		merge m:1 pid calendar_year using `everkia', nogen keep(1 3)

		foreach x in calendar_year everhfp everwia everkia {
			rename `x' `x'_`ytitle'
		}
		foreach x in everhfp {
			replace `x'_`ytitle' = 0 if `x'_`ytitle'==. & inrange(calendar_year_`ytitle',1992,2018)
			replace `x'_`ytitle' = . if calendar_year_`ytitle'<1992 | calendar_year_`ytitle'>2018
		}
		foreach x in everwia everkia {
			replace `x'_`ytitle' = 0 if `x'_`ytitle'==. & inrange(calendar_year_`ytitle',2001,2018)
			replace `x'_`ytitle' = . if calendar_year_`ytitle'<2001 | calendar_year_`ytitle'>2018
		}
	}

	drop calendar_year*
		
	egen any_hfp_0_19 = rowmax(everhfp_0 everhfp_1 everhfp_2 everhfp_3 everhfp_4 everhfp_5 everhfp_6 everhfp_7 everhfp_8 everhfp_9 everhfp_10 everhfp_11 everhfp_12 everhfp_13 everhfp_14 everhfp_15 everhfp_16 everhfp_17 everhfp_18 everhfp_19)
	egen any_wia_0_19 = rowmax(everwia_0 everwia_1 everwia_2 everwia_3 everwia_4 everwia_5 everwia_6 everwia_7 everwia_8 everwia_9 everwia_10 everwia_11 everwia_12 everwia_13 everwia_14 everwia_15 everwia_16 everwia_17 everwia_18 everwia_19)
	egen any_kia_0_19 = rowmax(everkia_0 everkia_1 everkia_2 everkia_3 everkia_4 everkia_5 everkia_6 everkia_7 everkia_8 everkia_9 everkia_10 everkia_11 everkia_12 everkia_13 everkia_14 everkia_15 everkia_16 everkia_17 everkia_18 everkia_19)

	foreach cmf in 11 12 13 19 25 42 68 88 92 94 {
		gen cmf_`cmf' = comcmf == `cmf'
	}
	gen othobc_cmf = conspmos!="" & (inlist(comcmf,11,12,13,19,25)!=1) & (inlist(comcmf,42,68,88,92,94)!=1)
	gen othcmf = conspmos!="" & (inlist(comcmf,11,12,13,19)!=1) & (inlist(comcmf,42,68,92,94)!=1)
	gen obs_mos = conspmos != ""
	tab conspmos if  inrange(firstyear,1990,2011) & inrange(firstafqt,12,49) & black == 1, sort
	tab conspmos if  inrange(firstyear,1990,2011) & inrange(firstafqt,12,49) & white == 1, sort
	tab conspmos if  inrange(firstyear,1990,2011) & inrange(firstafqt,31,68) & black == 1, sort
	tab conspmos if  inrange(firstyear,1990,2011) & inrange(firstafqt,31,68) & white == 1, sort
	
	*if doing top 3 from each
	gen mos92g = conspmos == "92g"
	gen mos11b = conspmos == "11b"
	gen mos92a = conspmos == "92a"
	gen mos91b = conspmos == "91b"
	gen mos88m = conspmos == "88m"
	gen mos42a = conspmos == "42a"
	gen mos31b = conspmos == "31b"

	*---- reduced form ----*
	eststo clear
	local round = 0
	foreach out in combat_tc noncombat_tc any_hfp_0_19 any_wia_0_19 any_kia_0_19 cmf_11 cmf_12 cmf_13 cmf_19 cmf_42 cmf_68 cmf_92 cmf_94 othcmf {
		local round = `round' + 1
		if "`out'" == "combat_tc" local lab "Combat Occupation Mean"
		if "`out'" == "combat_tc" local ys 1990

		if "`out'" == "noncombat_tc" local lab "Non-Combat Occupation"
		if "`out'" == "noncombat_tc" local ys 1990

		if "`out'" == "any_hfp_0_19" local lab "Ever Deployed (0-19)"
		if "`out'" == "any_hfp_0_19" local ys 1992

		if "`out'" == "any_wia_0_19" local lab "Ever WIA (0-19)"
		if "`out'" == "any_wia_0_19" local ys 2001

		if "`out'" == "any_kia_0_19" local lab "Ever KIA (0-19)"
		if "`out'" == "any_kia_0_19" local ys 2001

		if "`out'" == "cmf_11" local lab "Infantry"
		if "`out'" == "cmf_12" local lab "Corps of Engineers"
		if "`out'" == "cmf_13" local lab "Field Artillery"
		if "`out'" == "cmf_19" local lab "Armor"
		if "`out'" == "cmf_42" local lab "Adjutant General"
		if "`out'" == "cmf_68" local lab "Medical CMF"
		if "`out'" == "cmf_92" local lab "Quarter Master Corps"
		if "`out'" == "cmf_94" local lab "Maintenance"
		if "`out'" == "othcmf" local lab "Other CMF"

		gen inst = inst31  
		eststo `out'_b1: reg `out' inst k31 k31_2 instk31 instk31_2 i.quarterFE if inrange(firstyear,`ys',2011) & inrange(firstafqt,12,49) & black == 1, robust
		eststo `out'_w1: reg `out' inst k31 k31_2 instk31 instk31_2 i.quarterFE if inrange(firstyear,`ys',2011) & inrange(firstafqt,12,49) & white == 1, robust
		
		drop inst
		g inst = inst50
		eststo `out'_b2: reg `out' inst k50 k50_2 instk50 instk50_2 i.quarterFE if inrange(firstyear,`ys',2011) & inrange(firstafqt,31,68) & black == 1, robust
		eststo `out'_w2: reg `out' inst k50 k50_2 instk50 instk50_2 i.quarterFE if inrange(firstyear,`ys',2011) & inrange(firstafqt,31,68) & white == 1, robust
		drop inst
		
		noi esttab  `out'_b1  `out'_w1  `out'_b2  `out'_w2 ///
			using "${output}cumulative/occupation_table_panel`round'.tex",	///
			cells(b(fmt(3) star) se(par fmt(3)))  ///
			replace noobs nomtitle nodepvar nonumber plain fragment label ///	
			substitute(_ \_) style(tex) ///	
			mlabels(, none) collabels(, none) 	///			
			keep(inst, relax) ///
			varlabels(inst "`lab'") ///
			starlevels(* 0.1 ** 0.05 *** 0.01) varwidth(30) ///
			stats()
	}
			
	*---- 2SLS ----*
	eststo clear
	local round = 0
	foreach out in combat_tc noncombat_tc any_hfp_0_19 any_wia_0_19 any_kia_0_19 cmf_11 cmf_12 cmf_13 cmf_19 cmf_42 cmf_68 cmf_92 cmf_94 othcmf {
		local round = `round' + 1
		if "`out'" == "combat_tc" local lab "Combat Occupation"
		if "`out'" == "combat_tc" local ys 1990
		if "`out'" == "combat_tc" local endog obs_mos

		if "`out'" == "noncombat_tc" local lab "Non-Combat Occupation"
		if "`out'" == "noncombat_tc" local ys 1990
		if "`out'" == "noncombat_tc" local endog obs_mos

		if "`out'" == "any_hfp_0_19" local lab "Ever Deployed (0-19)"
		if "`out'" == "any_wia_0_19" local lab "Ever WIA (0-19)"
		if "`out'" == "any_kia_0_19" local lab "Ever KIA (0-19)"

		if inlist("`out'","any_hfp_0_19","any_wia_0_19","any_kia_0_19"){
			local ys 1992
			local qs 126
			local endog accessarmy
		}

		if "`out'" == "cmf_11" local lab "Infantry"
		if "`out'" == "cmf_11" local ys 1990
		if "`out'" == "cmf_11" local endog obs_mos

		if "`out'" == "cmf_12" local lab "Corps of Engineers"
		if "`out'" == "cmf_12" local ys 1990
		if "`out'" == "cmf_12" local endog obs_mos

		if "`out'" == "cmf_13" local lab "Field Artillery"
		if "`out'" == "cmf_13" local ys 1990
		if "`out'" == "cmf_13" local endog obs_mos

		if "`out'" == "cmf_19" local lab "Armor"
		if "`out'" == "cmf_19" local ys 1990
		if "`out'" == "cmf_19" local endog obs_mos

		if "`out'" == "cmf_42" local lab "Adjutant General"
		if "`out'" == "cmf_42" local ys 1990
		if "`out'" == "cmf_42" local endog obs_mos

		if "`out'" == "cmf_68" local lab "Medical CMF"
		if "`out'" == "cmf_68" local ys 1990
		if "`out'" == "cmf_68" local endog obs_mos

		if "`out'" == "cmf_92" local lab "Quarter Master Corps"
		if "`out'" == "cmf_92" local ys 1990
		if "`out'" == "cmf_92" local endog obs_mos

		if "`out'" == "cmf_94" local lab "Maintenance"
		if "`out'" == "cmf_94" local ys 1990
		if "`out'" == "cmf_94" local endog obs_mos

		if "`out'" == "othcmf" local lab "Other CMF"
		if "`out'" == "othcmf" local ys 1990
		if "`out'" == "othcmf" local endog obs_mos

		di "out=`out', endog=`endog', qs=`qs'"
		eststo `out'_b1: ivregress 2sls `out' k31 k31_2 instk31 instk31_2 i.quarterFE  (`endog'=inst31) if inrange(quarterFE,`qs',207) & inrange(firstafqt,12,49) & black == 1, robust
		eststo `out'_w1: ivregress 2sls `out' k31 k31_2 instk31 instk31_2 i.quarterFE  (`endog'=inst31) if inrange(quarterFE,`qs',207) & inrange(firstafqt,12,49) & white == 1, robust
		eststo `out'_b2: ivregress 2sls `out' k50 k50_2 instk50 instk50_2 i.quarterFE  (`endog'=inst50) if inrange(quarterFE,`qs',207) & inrange(firstafqt,31,68) & black == 1, robust
		eststo `out'_w2: ivregress 2sls `out' k50 k50_2 instk50 instk50_2 i.quarterFE  (`endog'=inst50) if inrange(quarterFE,`qs',207) & inrange(firstafqt,31,68) & white == 1, robust
		
		noi esttab  `out'_b1  `out'_w1  `out'_b2  `out'_w2 ///
			using "${output}cumulative/occupation_table_2sls_panel`round'_`endog'.tex",	///
			cells(b(fmt(3) star) se(par fmt(3)))  ///
			replace noobs nomtitle nodepvar nonumber plain fragment label ///	
			substitute(_ \_) style(tex) ///	
			mlabels(, none) collabels(, none) 	///			
			keep(`endog', relax) ///
			varlabels(`endog' "`lab'") ///
			starlevels(* 0.1 ** 0.05 *** 0.01) varwidth(30) ///
			stats(N, labels("\hline \\ Observations") fmt(0))
	}
}	//end if
else di "*---- skipping occupations ----*"

